# Replication code for Christopher Chambers-Ju & Mart Trasberg manuscript "Framing Labor Action" (Last updated Feb 15, 2025)
# Code for the Appendix A1 is in the companion R script

install.packages("foreign")
install.packages("haven")
install.packages("broom")
install.packages("dplyr")
install.packages("plyr")

# Loading data
library(foreign)
library(haven)

dta <- read_dta("C:/Users/L03538527/OneDrive - Instituto Tecnologico y de Estudios Superiores de Monterrey/Research/Education in Mexico/Drafts/Perspectives on Politics/Analysis/OA23 Base UT-TEC.dta")
#dta <- read_dta("OA23 Base UT-TEC.dta")
# Use "Import Dataset" in RStudio

# dta<-OA23_Base_UT_TEC
attach(dta)
objects()

### SUPPORT FOR UNIONS ###

# Relabeling Treatments

TratUT2<-as.factor(TratUT2)
levels(TratUT2)
library(dplyr)
library(plyr)
dta$Treatment<-revalue(TratUT2, c("1"="Control","2"="Strike","3"="Pro_strike",
                                    "4"="Anti_strike","5"="Double_frame"))

# Relabeling union support DV
# "Pensando en los sindicatos en general, ¿le gustaría personalmente que los sindicatos en México tuvieran más influencia, la misma influencia o menos influencia que la que tienen hoy?

UT2_1<-as.numeric(UT2_1)
dta$unions <- dta$UT2_1

# Recode the new variable, so higher values are stronger support for unions
dta$unions <- ifelse(dta$unions == 1, 5,
                      ifelse(dta$unions == 2, 4,
                             ifelse(dta$unions == 3, 3,
                                    ifelse(dta$unions == 4, 2,
                                           ifelse(dta$unions == 5, 1, NA)))))

# New variable checks with the old variable
table(dta$unions,dta$UT2_1)

# Descriptive exploration, mean of U2_1 for each category
mean_unions_Tr_2<-dta   %>% 
  group_by(dta$Treatment) %>%
  dplyr::summarize(Mean = mean(unions, na.rm=TRUE))  

mean_unions_Tr_2
# Mean Control for unions = 3.25

# Descriptive exploration, calculating the mean and confidence intervals for 
# each treatment on "union support"

summary_expertiment_union <- dta %>%
  group_by(Treatment) %>%
  dplyr::summarise(mean_outcome = mean(unions, na.rm = TRUE),
            ci_lower = mean_outcome - qt(0.95, length(unions)-1) * sd(unions,na.rm = TRUE) / sqrt(length(unions)),
            ci_upper = mean_outcome + qt(0.95, length(unions)-1) * sd(unions,na.rm = TRUE) / sqrt(length(unions)))

summary_expertiment_union


### FIGURE 1 ####

library(MASS)
library(broom)
library(dplyr)
library(ggplot2)

# Fit the ordered logistic regression model and tidy the output with 95% confidence intervals
tidy_reg_conf_int <- dta %>% 
  polr(as.factor(unions) ~ Treatment, data = ., Hess = TRUE) %>%  # Fit ordered logistic model
  tidy(conf.int = TRUE, conf.level = 0.95)  # Specify 95% confidence intervals

# Display the tidy data with confidence intervals
tidy_reg_conf_int

# Rename variables and customize their order
tidy_reg_conf_int <- tidy_reg_conf_int %>%
  mutate(term = recode(term,
                       "TreatmentControl" = "Control",
                       "TreatmentStrike" = "Strike",
                       "TreatmentPro_strike" = "Pro-Strike",
                       "TreatmentAnti_strike" = "Anti-Strike",
                       "TreatmentDouble_frame" = "Double Frame")) %>%
  # Convert term to a factor and customize the order
  mutate(term = factor(term, levels = c("Double Frame", "Anti-Strike", "Pro-Strike", "Strike", "Control")))

# Create the coefficient plot with error bars for 95% confidence intervals
tidy_reg_conf_int %>%
  filter(term != "(Intercept)") %>%  # Exclude the intercept term
  ggplot(aes(estimate, term)) +  # Different colors for different confidence levels
  geom_point(position = position_dodge(width = 0.5)) +  # Slightly dodge points to avoid overlap
  geom_errorbarh(aes(xmin = conf.low, xmax = conf.high), 
                 position = position_dodge(width = 0.5), 
                 height = 0.2) +  # Add horizontal error bars with slight dodge
  geom_vline(xintercept = 0, lty = 2) +  # Add a dotted vertical line at zero
  labs(
    x = "Estimate of treatment effect on union support",
    y = NULL,
    title = "Coefficient plot with 95% confidence intervals",
    color = "Confidence Level"
  ) +
  theme_minimal()  # Apply a minimal theme for a cleaner look



### FIGURE 2 ###

dta$Treatment <- relevel(dta$Treatment, ref = "Pro_strike")
summary(lm(dta$unions ~ dta$Treatment))
# Create a new variable with the desired reference category
#dta$new_Treatment <- factor(dta$Treatment, levels = c("desired_reference_category", "other_levels"))
# Run the linear model with the new variable
#summary(lm(dta$unions ~ dta$new_Treatment))

# Fit the ordered logistic regression model and tidy the output with 95% confidence intervals
tidy_reg_conf_int <- dta %>% 
  polr(as.factor(unions) ~ Treatment, data = ., Hess = TRUE) %>%  # Fit ordered logistic model
  tidy(conf.int = TRUE, conf.level = 0.95)  # Specify 95% confidence intervals

# Display the tidy data with confidence intervals
tidy_reg_conf_int

# Rename variables and customize their order
tidy_reg_conf_int <- tidy_reg_conf_int %>%
  mutate(term = recode(term,
                       "TreatmentControl" = "Control",
                       "TreatmentStrike" = "Strike",
                       "TreatmentPro_strike" = "Pro-Strike",
                       "TreatmentAnti_strike" = "Anti-Strike",
                       "TreatmentDouble_frame" = "Double Frame")) %>%
  # Convert term to a factor and customize the order
  mutate(term = factor(term, levels = c("Double Frame", "Anti-Strike", "Pro-Strike", "Strike", "Control")))

# Create the coefficient plot with error bars for 95% confidence intervals
tidy_reg_conf_int %>%
  filter(term != "(Intercept)") %>%  # Exclude the intercept term
  ggplot(aes(estimate, term)) +  # Different colors for different confidence levels
  geom_point(position = position_dodge(width = 0.5)) +  # Slightly dodge points to avoid overlap
  geom_errorbarh(aes(xmin = conf.low, xmax = conf.high), 
                 position = position_dodge(width = 0.5), 
                 height = 0.2) +  # Add horizontal error bars with slight dodge
  geom_vline(xintercept = 0, lty = 2) +  # Add a dotted vertical line at zero
  labs(
    x = "Estimate of treatment effect on union support",
    y = NULL,
    title = "Coefficient plot with 95% confidence intervals",
    color = "Confidence Level"
  ) +
  theme_minimal()  # Apply a minimal theme for a cleaner look




## SUPPORT for STRIKING TEACHERS
 
# Recoding support for striking teachers

table(dta$UT2_2) # 1 = totally oppose, 7 = totally support, 8 = NA
dta$Support_teachers<-ifelse(UT2_2 == 8, NA, UT2_2)
table(dta$Support_teachers)
class(dta$Support_teachers)

# Calculating the mean and confidence intervals for each treatment on "support for striking teachers"

summary_expertiment_teachers <- dta %>%
  group_by(Treatment) %>%
  dplyr::summarise(mean_outcome = mean(Support_teachers, na.rm = TRUE),
                   ci_lower = mean_outcome - qt(0.975, length(Support_teachers)-1) * sd(Support_teachers,na.rm = TRUE) / sqrt(length(Support_teachers)),
                   ci_upper = mean_outcome + qt(0.975, length(Support_teachers)-1) * sd(Support_teachers,na.rm = TRUE) / sqrt(length(Support_teachers)))

summary_expertiment_teachers

### FIGURE 3 ### 

## Control as reference category
dta$Treatment <- relevel(dta$Treatment, ref = "Control")
summary(lm(dta$Support_teachers ~ dta$Treatment))

# Load required libraries
library(MASS)
library(broom)
library(dplyr)
library(ggplot2)

# Fit the ordered logistic regression model and tidy the output with 95% confidence intervals
tidy_reg_conf_int <- dta %>% 
  polr(as.factor(Support_teachers) ~ Treatment, data = ., Hess = TRUE) %>%  # Fit ordered logistic model
  tidy(conf.int = TRUE, conf.level = 0.95)  # Specify 95% confidence intervals

# Rename variables and customize their order
tidy_reg_conf_int <- tidy_reg_conf_int %>%
  mutate(term = recode(term,
                       "TreatmentControl" = "Control",
                       "TreatmentStrike" = "Strike",
                       "TreatmentPro_strike" = "Pro-Strike",
                       "TreatmentAnti_strike" = "Anti-Strike",
                       "TreatmentDouble_frame" = "Double Frame")) %>%
  # Convert term to a factor and customize the order
  mutate(term = factor(term, levels = c("Double Frame", "Anti-Strike", "Pro-Strike", "Strike", "Control")))

# Create the coefficient plot with error bars for 95% confidence intervals
tidy_reg_conf_int %>%
  filter(term != "(Intercept)") %>%  # Exclude the intercept term
  ggplot(aes(estimate, term, color)) +  # Different colors for different confidence levels
  geom_point(position = position_dodge(width = 0.5)) +  # Slightly dodge points to avoid overlap
  geom_errorbarh(aes(xmin = conf.low, xmax = conf.high), 
                 position = position_dodge(width = 0.5), 
                 height = 0.2) +  # Add horizontal error bars with slight dodge
  geom_vline(xintercept = 0, lty = 2) +  # Add a dotted vertical line at zero
  labs(
    x = "Estimate of treatment effect on union support",
    y = NULL,
    title = "Coefficient plot with 95% confidence intervals",
    color = "Confidence Level"
  ) +
  theme_minimal()  # Apply a minimal theme for a cleaner look


### FIGURE 4 ###

##  Show Pro-strike as reference category
dta$Treatment <- relevel(dta$Treatment, ref = "Pro_strike")
summary(lm(dta$Support_teachers ~ dta$Treatment))

# Fit the ordered logistic regression model and tidy the output with 95% confidence intervals
tidy_reg_conf_int <- dta %>% 
  polr(as.factor(Support_teachers) ~ Treatment, data = ., Hess = TRUE) %>%  # Fit ordered logistic model
  tidy(conf.int = TRUE, conf.level = 0.95)  # Specify 95% confidence intervals

# Rename variables and customize their order
tidy_reg_conf_int <- tidy_reg_conf_int %>%
  mutate(term = recode(term,
                       "TreatmentControl" = "Control",
                       "TreatmentStrike" = "Strike",
                       "TreatmentPro_strike" = "Pro-Strike",
                       "TreatmentAnti_strike" = "Anti-Strike",
                       "TreatmentDouble_frame" = "Double Frame")) %>%
  # Convert term to a factor and customize the order
  mutate(term = factor(term, levels = c("Double Frame", "Anti-Strike", "Pro-Strike", "Strike", "Control")))

# Create the coefficient plot with error bars for 95% confidence intervals
tidy_reg_conf_int %>%
  filter(term != "(Intercept)") %>%  # Exclude the intercept term
  ggplot(aes(estimate, term, color)) +  # Different colors for different confidence levels
  geom_point(position = position_dodge(width = 0.5)) +  # Slightly dodge points to avoid overlap
  geom_errorbarh(aes(xmin = conf.low, xmax = conf.high), 
                 position = position_dodge(width = 0.5), 
                 height = 0.2) +  # Add horizontal error bars with slight dodge
  geom_vline(xintercept = 0, lty = 2) +  # Add a dotted vertical line at zero
  labs(
    x = "Estimate of treatment effect on union support",
    y = NULL,
    title = "Coefficient plot with 95% confidence intervals",
    color = "Confidence Level"
  ) +
  theme_minimal()  # Apply a minimal theme for a cleaner look


